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ABSTRACT 

Increasing  the  nitramine  content  of  solid  rocket  propellants  increases 
the  overall  performance  of  thp  system  as  well  as  the  sensitivity  to  Shock  to 
Detonation  Transition  (SDT)  l  Deflagration  to  Detonation  Transition  (DDT). 
This  report  deals  primarily  with  the  analysis  and  numerical  modeling  of  a 
combined  SDT/DDT  event.  The  results  show  that  in  some  instances  a  zone  of 
burning  granulated  propellant,  confined  and  adjacent  to  a  zone  of  cast 
propellant,  can  provide  a  rapid  enough  pressure-rise  rate  to  shock  initiate 
the  cast  material.  This  type  of  detonation  hazard  scenario  is  a  real 
possibility  in  any  high-energy  rocket  motor  environment. 

The  modeling  study  also  indicates  areas  where  important  assumptions  need 
to  be  further  researched.  These  include:  (a)  relations  for  dynamic  (tran¬ 
sient)  collapse  of  the  voids  or  pores;  (b)  relations  for  setting  the  volume 
percent  of  hot  spots  based  on  initial  porosity;  (c)  the  evaluation  and 
expression  for  the  chemical  rate  of  decomposition  of  the  reactive,  shocked 
material;  and  (d)  the  assessment  of  two-phase  mixture  equilibrium. 

The  predicted  run-to  detonation  distance  as  a  function  of  porosity  for 
HMX  explosive  compares  favorably  with  limited  shock  initiation  experiments. 
There  is  no  data  available  to  check  whether  the  predictions  of  ramp-wave  com¬ 
pressions  (where  rise  times  exceed  several  microseconds)  presented  here  are 
valid. 
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CHAPTER  I 

TRANSITION  TO  OETONATION:  THE  MODEL 


Introduction 

In  1984  a  model  was  developed  yielding  preliminary  results  that  predicted 
the  processes  of  transition  to  detonation  by  ramp-wave  compression  of  porous 
explosives  [1].  This  report  highlights  the  major  features  of  this  model  and 
presents  additional  analysis  which  we  believe  indicates  that  the  physical  and 
chemical  processes  in  the  model  are  properly  described. 

°redicting  Deflagration  to  Detonation  Transition  (DDT)  depends  on  being 
able  to  prescribe  the  shock  loading  conditions,  i.e.  the  rate  r.t  which  an  ex¬ 
plosive  is  compressed.  An  actual  shock  loading  condition  represents  sub¬ 
microsecond  rise  times,  while  ramp-wave  loadings  could  compress  the  solid  in 
relatively  long  times,  exceeding  tens  of  microseconds. 

The  motivation  of  our  work  is  to  assess  the  hazard  of  DDT  in  high  energy 
solid  propellant  rocket  motors.  Data  and  analysis  indicate  that  sufficient 
porosity  in  the  propellant  is  necessary.  Thus  in  order  for  a  DDT  to  occur  in 
a  full-scale  motor,  the  grain  must  first  be  damaged.  This  has  been  demon¬ 
strated  in  the  laboratory  many  times  when  small-scale  (20  to  40  cm)  packed 
beds  of  granulated  propellant  undergo  the  transition  from  deflagration  to 
detonation.  Combustion  of  the  high  surface-to-volume  ratio  micron-sized  frag¬ 
ments  provides  the  rapid  pressurization  rate  necessary  to  shock  initiate  the 
remaining  fragments.  By  increasing  the  nitramine  content  of  the  propellant 
mixture,  the  propensity  to  detonate  is  increased  in  two  ways:  the  decompos¬ 
ition  rate  of  the  propellant  increases  and  the  shock  sensitivity  increases. 
Although  little  has  been  done  to  experimentally  verify  it,  one  can  also  assume 
that  the  particle  morphology  (size,  shape,  and  surface  roughness)  greatly  in¬ 
fluences  the  decomposition  rate  and  thus  the  build  up  to  detonation. 

We  have  proposed  two  possible  methods  whereby  DDT  can  occur  in  a  solid 
rocket  motor.  In  the  first  case  a  length,  L,  of  the  motor  grain  is  completely 
fragmented.  This  case  was  the  topic  of  several  research  papers  [ 2-4 1  by  our 
group  and  will  be  discussed  here  only  to  provide  a  comparison  with  the  second 
case  proposed.  In  Case  I,  the  bed  is  of  sufficient  length  to  allow  the  trans¬ 
ition  to  detonation  to  occur  within  the  granulated  material,  that  is,  L  >  L$gy 
where  L^qj  is  the  detonation  run-up  length  measured  from  the  ignition 
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source.  In  Case  2,  L  <  L^qj  and  the  detonation  occurs  in  the  upstream  cast 
material.  However  in  both  cases,  the  confined  burning  of  the  propellant  frag¬ 
ments  is  what  drives  the  deflagration  wave  to  a  high-order  detonation.  The 
reader  should  consult  Ref.  [5]  for  a  general  overview  of  DDT  in  different  pro¬ 
pellants  and  explosives  and  Ref.  [6]  for  typical  experimental  results. 

In  Case  1,  the  length  of  the  packed  bed  is  longer  than  the  critical 
length  necessary  for  accelerated  convective  combustion  to  occur,  and  sub¬ 
sequent  detonation  is  inevitable.  Figure  1  illustrates  the  Case  1  bed  con¬ 
figuration.  The  second  DDT  scenario  involves  a  region  of  burning  granulated 
propellant  providing  the  impetus  to  shock  initiate  an  adjacent  region  of  the 
case  propellant.  This  is  illustrated  in  Fig.  2.  The  cast  material  (Zone  1) 
can  contain  "blind"  pores  and  is  assumed  to  be  impermeable  to  the  flow  of  hot 
gases  from  the  zone  of  fragmented  material  (Zone  2).  This  implies  that,  un¬ 
like  the  first  DDT  scenario  discussed,  only  stress  waves  can  be  transmitted 
across  the  Zone  1/Zone  2  interface.  Figure  3  is  a  schematic  representation  of 
the  sequence  of  events  leading  to  this  type  of  DDT.  Superimposed  on  each  sec¬ 
tion  of  the  figure  is  a  solid  line  representing  the  local  gas  porosity  (volume 
of  gas/total  volume)  as  a  function  of  x,  the  bed  location.  A  value  of  <t>  equal 
to  unity  represents  a  zone  of  all  gas  while  <j>  equal  to  zero  indicates  a  void¬ 
less  solid. 

Part  A  shows  a  zone  of  burning  granular  propellant  adjacent  to  a  zone  of 
porous,  cast  propellant.  Here,  the  heavy  black  dots  are  representative  of 
microvoids  in  the  cast  material.  Illustrated  in  Part  B  is  the  collapse  of  the 
pores,  a  result  of  the  stress  load  transmitted  across  the  granulated  bed/cast 
propellant  interface.  Parts  C  and  D  show  the  length  of  the  pore  collapse  zone 
to  increase  with  time  as  the  lead  compression  waves  travel  farther  into  the 
propellant  bed.  The  finite  compression  waves  coalesce  into  a  shock  front 
which  subsequently  shock  initiates  the  cast  propellant  at  a  location  ahead 
of  the  interface.  From  this  location  a  detonation  wave  propagates  through  the 
porous  material  while  a  retonation  wave  propagates  back  through  the  compressed 
material  (Part  E). 

In  related  research,  Setchell  [71  studied  ramp-wave  initiation  mechanisms 
using  waves  with  rise-times  of  0.3  to  0.8  ys.  In  Ref.  (7]  the  ramp-waves  were 
developed  by  propagation  of  an  impact-generated  shock  wave  through  a  unique 
material  which  has  the  property  of  spreading  out  the  shock  front  into  a  ramp- 


Figure  1  Enlargement  of  granulated  bed  showing  processes 
that  need  to  be  modeled  for  DSDT. 
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wave.  In  the  work  presented  here,  the  rise-times  were  much  longer  than  those 
cited  in  [91,  typically  of  the  order  of  tens  of  microseconds.  These  slower 
rise-times  are  typical  pressure-time  histories  generated  within  a  deflagrating 
bed  of  granular  HMX. 

The  purpose  of  this  work  is  to  model  the  key  elements  of  the  six-part 
scenario  shown  in  Fig.  3.  A  one-dimensional  hydrodynamic  Lagrangian  finite 
difference  technique  is  used  to  numerically  solve  the  conservation  equa¬ 
tions.  A  static  pore  collapse  model  [8]  which  demarcates  three  regimes  of 
deformation— elastic,  elastic-plastic,  and  plastic— is  utilized  to  determine 
the  rate  of  compaction  and  the  development  of  the  solid  plug.  Chapter  3 
investigates  two  dynamic  pore  collapse  models.  Furthermore,  a  hot  spot  modal 
[9]  is  incorporated  in  the  code  to  define  the  sensitivity  to  reaction.  By 
introducing  reactive  chemistry  to  the  code,  a  strong  effort  is  made  to  model 
the  detonation  and  retonation  waves  which  are  initiated  by  a  shock  wave 
generated  from  ramp-wave  inputs  with  rise  times  on  the  order  of  tens  of 
microseconds.  It  was  shown  in  Ref.  [11  that  rise  times,  t*,  of  this  order  are 
typical  for  burning,  granulated  beds  which  have  lengths  less  than  their 
critical  detonation  run-up  length. 

The  Model 

The  scope  of  the  work  presented  in  this  report  is  to  model  DDT  in  the 
granulated  bed/cast  propellant  configuration  shown  in  Fig.  2.  It  is  assumed 
that  the  granulated  bed  is  not  long  enough  (L  <  L$gy  to  undergo  DDT  in  the 
usual  sense  (Case  1)  but,  by  convective  burning,  can  provide  the  driving  force 
necessary  to  shock  initiate  the  upstream  cast  propellant.  For  this  analysis 
the  cast  propellant  is  assumed  to  have  "blind"  pores  and  to  be  impermeable  to 
the  flow  of  product  gases  generated  in  the  granular  bed.  This  implies  that 
only  stress  waves  can  be  transmitted  across  the  granular  bed/cast  propellant 
interface.  The  rate  at  which  the  reacting  granulated  bed  stresses  the  cast 
propellant  is  a  function  of  many  of  the  granular  bed  parameters  including  par¬ 
ticle  size,  loading  density,  solid  density,  and  reaction  rate.  In  the 
analysis  to  follow,  the  boundary  pressure-rise  rate  was  determined  by 
executing  the  DDT  code  discussed  in  Ref.  [31  for  granular  beds  with  lengths, 
L,  less  than  the  detonation  run-up  length  and  recording  the  pressure-rise  rate 
at  the  end  opposite  the  Igniter.  The  pressurization  rate  in  the  granular  bed, 
dP*/dt*,  is  strongly  dependent  on  the  size  of  the  particles  being  consumed. 
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The  larger  the  particle  being  consumed,  the  slower  the  pressurization  rate  on 
the  interface.  Thus,  from  this  type  of  analysis  we  have  been  able  to 
determine  the  P(t)  boundary  condition  for  the  stress  wave  analysis  of  the 
upstream  cast  propellant. 

The  P(t)  functions  obtained  from  this  modeling  effort  have  been 
linearized  for  this  work  and  are  expressed  as 

P(t)  -  (P*-P0)(t/t*)+  Pc  t  <  t*  (1.1a) 

and 

P(t)  =  P*  t  <  t*  ( 1 . Id ) 

where  P*  is  the  maximum  stress  in  the  bed.  The  parameter  t*  is  a  characteris¬ 
tic  rise-time  for  the  ramp-wave  input  function. 

The  Lagrangian  or  material  form  of  the  governing  equations  are  in¬ 
corporated  in  the  hydrodynamic  analysis  of  a  continuous  material  with  a  moving 
boundary.  In  the  problem  addressed,  the  moving  boundary  is  a  result  of  the 
applied  stress  load  from  the  burning  granulated  bed.  The  one-dimensional  un¬ 
steady  form  of  the  conservation  of  mass,  momentum,  and  energy  equations  are 
expressed  for  the  two-phase  mechanical  mixture  as 


vt  = 

;  vux’ 

(1.2) 

ut  = 

1 

< 

"O 

X 

* 

(1.3) 

and 

et  = 

!  -Pvt  +  Q. 

(1.4) 

In 

the 

above  expressions,  v  represents  the 

specific 

volume  of 

the 

mixture; 

u. 

the  particle  velocity;  e. 

the 

specific 

internal 

energy;  P, 

the 

total  stress 

;  and  Q,  the  heat  added 

per 

unit  mass  per  unit  time. 

The 

subscripts  x  and  t  indicate  partial  derivatives  with  respect  to  the  Lagrangian 
spatial  coordinate  and  time,  respectively.  In  addition  to  the  conservation 
equations,  a  material  equation  of  state,  P  =  P(v,e),  and  appropriate 
constitutive  relations  provide  for  mathematical  closure. 

For  the  solid  material  the  equation  of  state  P$(vs,es)  and  the  caloric 
equation  e$(vs,Ts)  are  expressed  in  terms  of  a  Helmholtz  free  energy  function 
[10|,  i|i(vs,T),  and  its  thermodynamic  derivatives  through  the  Second  Law  of 
Thermodynamics  reciprocity  relations 


a 


Ps  =  -  (H/3Vs)t 
es  =  <i  -  T$  (34</3T $)v. 


(1.5) 


(1.6) 


With  the  assumption  that  the  Gruneisen  coefficient  is  constant,  the  Helmholtz 
free  energy  function  takes  the  form 

*  -  J(vs)  *  rC„sln<v„/vs)<YTo>  +  Cvs[Tsln(T0/Ts)  +  T^l  (1.7) 

where  r  is  the  Gruneisen  coefficient  defined  by  the  thermodynamic  derivative 


r(v)  -  -  v  (£)„ 


(1.8) 


and  Cvs  is  the  specific  heat  at  constant  volume  of  the  solid  phase 

C  =  (—) 
vs  4tV 


(1.9) 


The  term  J(v$)  in  Eq.  (1.7)  is  a  nonlinear  volume-dependent  function 
determined  from  shock  Hugoniot  experiments  (10]. 

With  the  introduction  of  product  gases  into  the  system,  a  constitutive 
law  for  the  gas  phase  must  also  be  provided.  A  non-ideal  covolume  equation  of 
state  is  utilized 


Pg  -  -W1  +  6“g> 


(1.10) 


where  R  is  the  gas  constant  and  6  is  a  covolume  correction  term.  The  value 
of  s  is  determined  from  the  values  for  pressure,  temperature,  and  density  at 
the  CJ  state  predicted  by  the  TIGER  code  [11].  Table  1  gives  a  listing  of 
those  values  for  several  loading  densities  of  HMX.  In  accordance  with  the  re¬ 
ciprocity  relations  defined  earlier  (Eqs.  (5)  and  (6)),  the  caloric  equation 
of  state  for  the  gas  phase  is 


deg  =  CvgdTg 


(1.11) 


where  Cyg  is  the  specific  heat  at  constant  volume  of  the  product  gases. 
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Table  1  CJ  Parameters  (P,T,v,D,  from  Tiger) 


2o_ 

PCJ  (GPa) 

TCJ  <K) 

VCJ  (cc/9) 

D  (mm/us) 

6  (cc/g) 

1 

36.57 

3840 

0.4067 

9.20 

4.86 

1.06 

32.25 

3973 

0.4260 

8.77 

4.50 

1.12 

28.50 

4083 

0.4471 

8.36 

4.22 

1.19 

25.20 

4171 

0.4707 

7.98 

4.00 

1.27 

22.22 

4280 

0.5276 

7.33 

3.74 

1.46 

17.26 

4304 

0.5629 

7.05 

3.68 
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The  governing  equations  (Eqs.  (2)  through  (4))  are  expressed  In  terms  of 
the  thermodynamic  properties  (P,v,e)  and  dynamic  variable  u  of  the  two-phase 
mixture.  The  Individual  phase  volumes  sum  to  equal  the  total  Vy  =  +  \l^ 

where  =  v^m^  (i  =  g,s).  Likewise  for  the  mixture  total  energy 

ET  3  EG  +  ES* 

We  define  the  material  porosity  (a)  as  the  ratio  of  total  volume  to 
volume  occupied  by  the  solid  phase 


a  -  VT/VS 


(1.12) 


and  incorporate  a  three-step  static  material  collapse  model  [81.  Current  work 
by  our  research  group  involves  an  investigation  of  the  more  realistic  dynamic 
material  response  models. 

A  first  order  Arrhenius  burn  model  is  used  to  describe  the  chemical  de¬ 
composition  of  the  reactive  material 


^  =  -ZW  exp(-E*/RT*)  . 
dt 


(1.13) 


In  the  above  equation,  W  is  the  mass  fraction  of  unreacted  explosive, 
W  =  ms/m$0;  Z  is  the  frequency  factor,  E*  the  activation  energy,  R  the  uni¬ 
versal  gas  constant,  and  T*  the  characteristic  burn  temperature.  During  com¬ 
pression  of  the  porous  bed,  T*  represents  a  "hot  spot"  temperature  due  to  ir¬ 
reversible  heating  at  pore  sites.  The  localized  hot  spot  temperature  is  dif¬ 
ferent  from  the  bulk  shock  temperature  and  is  calculated  from  an  energy  part¬ 
ition  model  [9[ .  The  underlying  assumption  in  this  hot-spot  model  is  that  the 
shocked  porous  material  is  at  one  of  two  possible  temperatures,  a  bulk  shock 
temperature,  T,  or  a  hot  spot  temperature,  TH.  The  energy  deposited  by  the 
shock  wave  is  equated  on  a  mass  fraction  basis  to  the  sum  of  the  reversible 
work  done  In  isentroplcally  compressing  the  bulk  of  the  material  plus  the  ir¬ 
reversible  heating  of  localized  hot  spots  [91 


<Vto-V>  "  WHe<v’TH>  + 


(1.14) 


In  Eq.  (1.14),  the  left-hand  side  represents  the  total  energy  deposited  in  the 
material  by  the  shock  of  strength  P.  The  term  e^(P)  represents  the  energy  re- 


quired  to  isentropical ly  compress  the  bulk  of  the  material  to  the  final  shock 
pressure  and  the  remaining  energy  term,  e(v,TH),  is  the  energy  available  to 
irreversibly  heat  the  hot  spots.  This  model  assumes  the  mass  fraction  of  the 
hot  spots,  W^j ,  to  be  equal  to  the  preshock  volume  fraction  of  pores. 

WH  (vto/vso_1)  =  *o  -  1  f1*15) 

Here,  the  subscript  "to"  represents  the  initial  porous  state  and  the  subscript 
"so"  refers  to  the  homogeneous  initial  state.  This  particular  macroscopic  hot 
spot  model  was  chosen  since  it  is  more  applicable  to  the  problem  at  hand.  The 
reader  should  consult  the  literature  for  one  of  the  more  detailed  microscopic 
approaches  to  the  theory  of  hot  spots. 

Equations  (1.2)  through  (1.15)  define  the  fluid  motion  and  thermodynamic 
state  of  continuous,  reactive,  porous  media.  The  system  of  equations  is 
solved  by  using  a  finite  difference  numerical  technique.  At  t  =  0,  the  bed  of 
porous  propellant  is  discretized  into  J  cells  labeled  from  left  to  right  as 

j  =  1,2 . J.  The  thermodynamic  properties  (pressure,  temperature,  internal 

energy,  and  specific  volume)  are  assumed  to  be  constant  over  the  width  of  each 
cell.  At  the  boundaries  between  the  cells,  values  for  particle  velocity  and 
spatial  location  are  assumed  known  at  t  =  0.  The  reader  is  referred  to 
Ref.  (12]  for  a  listing  of  the  finite  difference  approximations  to  the  govern¬ 
ing  differential  equations  and  constitutive  relations. 

Results 

As  stated  earlier,  the  purpose  of  this  work  is  to  predict  the  transient 
events  leading  to  detonation  in  a  bed  of  porous,  cast  propellant  being  stress 
loaded  at  one  end  by  a  burning,  granulated  bed  of  the  same  material.  In  our 
preliminary  calculations,  we  found  that  rise-times  of  10  ys  <  t*  <  30  ys  were 
typical  of  most  of  the  granulated  beds  studied.  Note  that  these  are  greater 
than  those  cited  earlier  [7], 

Several  cases  with  various  maximum  input  stresses,  characteristic  rise 
times,  and  initial  porosities  were  chosen  to  test  the  trends  of  the  computer 
model.  The  first  case  considered  has  an  initial  porosity  of  aQ  =  1.1176  and 
is  initiated  by  a  "ramp"  wave  with  a  maximum  input  pressure  of  P*  =  2  GPa  and 
characteristic  rise  time  of  t*  =  10  ysec.  Figures  4  through  6  illustrate  the 


12 


distribution  history  of  three  key  parameters,  namely  P,  a,  and  (1-W).  For 
this  particular  case,  three  specific  segments  in  time  were  chosen  for  il¬ 
lustration  purposes.  The  first  time  shown  in  all  the  three  separate  figures 
is  t  =  10  psec.  By  viewing  Fig.  4a,  one  can  see  fhat  the  compression  wave  has 
propagated  into  the  bed  to  a  distance  of  1.5  cm.  Notice  even  though  the  left 
boundary  has  reached  the  maximum  input  stress  shown  approximately  at  a  loca¬ 
tion  of  2  mm  into  the  bed,  the  shock  front  has  not  fully  developed.  Figure  5a 
illustrates  the  closure  of  the  voids  as  a  result  of  the  compression  wave.  One 
can  see  from  Fig.  6a  (gas-phase  mass  fraction)  that  the  reaction  has  yet  to 
commence. 

As  time  progresses  to  12  psec,  Fig.  6b  shows  that  the  reaction  has  com¬ 
menced.  This  occurrence  can  also  be  viewed  by  an  increase  in  porosity  (a),  as 
displayed  in  Fig.  5b  at  that  instant.  Induced  by  the  initiation  of  propellant 
decomposition,  an  increase  in  strength  of  the  compression  front  occurs,  as  is 
evident  in  Fig.  4b.  With  the  increase  in  the  degree  of  reaction,  the  chemical 
energy  further  contributes  to  strengthen  the  compression  front. 

At  t  =  13  usee,  the  decomposition  of  the  propellant  is  complete  in  a 
small  region  of  the  bed  near  x  =  2  cm,  as  shown  in  Fig.  6c.  The  porosity  dis¬ 
tribution,  displayed  in  Fig.  5c  (at  13  psec),  shows  compression  ahead  of  the 
decomposition  region.  The  net  effect  of  all  of  this  is  a  right  movi-.;j  detona¬ 
tion  wave  and  a  left  (or  rearward)  moving  retonation  wave,  as  illustrated  in 
Fig.  4c.  Note  that  the  retonation  wave  has  a  higher  peak  pressure  than  the 
detonation  wave.  This  is  due  to  the  higher  density  of  the  material  the  wave 
is  propagating  through.  The  steady  detonation  wave  continues  to  propagate 
even  further  into  the  bed  at  a  predicted  CJ  pressure  of  25  GPa  and  a  corres¬ 
ponding  CJ  temperature  of  3923  K. 

In  addition  to  the  figures  just  presented,  the  locus  of  the  stress  and 
reaction  fronts  are  shown  in  Fig.  7.  The  dashed  line  depicts  the  stress  wave 
propagation  into  the  bed,  initiating  the  detonation  at 
t  =  12.41  psec,  x  =  2.05  cm.  The  solid  lines  represent  the  locus  of  the  right 
and  left  traveling  detonation  fronts.  Figure  7  also  shows  the  termination 
point  of  the  retonation  wave.  A  change  in  velocity  is  apparent  in  the  left 
traveling  wave  at  the  termination  point.  The  slopes  of  the  solid  lines  cor¬ 
respond  to  a  detonation  and  retonation  velocity  of  D  =  8.8  mm/usec  and 
R  =  9.1  mn/psec,  respectively. 
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Figure  4  (a)  Pressure-distance  profiles  for  porous  HMX. 

Other  inputs  are  P*  =  2  GPA,  t*  =  10  psec. 
At  t  =  10  ysec,  no  reaction  is  noticeable, 
(b)  At  t  -  12  psec,  reaction  has  commenced. 


Pressure-distance  profile  at  t  =  13  sec  showing 
full  reaction  with  PCJ  =  25  GPA,  D  =  8.8  mm/psec 
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Porosity-distance  profile  at  t  =  13  wsec 
showing  detonated  zone  (1.6  <  x  <  2.6  cm) 
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reaction  is  noticeable. 
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F'gure  7  Predicted  locus  of  the  compression  front  with  transition 
to  both  detonation  and  a  rearward  moving  retonation 
(aQ  =  1.1176,  P*  =  2  GPa,  t*  -  10  wsec). 
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As  stated  earlier  in  the  text,  the  granular  bed  (Zone  2,  Fig.  2)  provides 
the  driving  force  to  shock  initiate  the  cast  material  (Zone  1,  Fig.  2).  In 
addition,  the  pressure-rise  rate  at  the  Zone  I/Zone  2  interface  is  a  function 
of  several  key  parameters  in  the  granulated  bed.  Figure  8  illustrates  the  in¬ 
fluence  of  pressure-rise  rate  on  the  run-up  distance  to  detonation,  L^gy. 
Data  are  presented  for  four  different  values  of  t*  (1,  3,  5,  and  10  ys, 
respectively).  The  slower  rise  time  of  10  ys  corresponds  to  a  granular  bed 
with  d0  =  500  ym  particles  driving  the  shock  build-up  process.  Previous 
calculations  show  that  HMX  particles  of  this  size  would  require  a  rather  long 
bed  length  in  order  to  exhibit  DDT  in  the  Case  1  situation.  In  general,  one 
would  expect  a  greater  run-up  distance  to  correspond  to  a  longer  characteris¬ 
tic  rise  time  for  a  specific  input  pressure,  since  the  compression  front  will 
reach  the  critical  pressure  needed  to  initiate  detonation  later.  Indeed, 
Fig.  8  demonstrates  the  occurrence  of  a  longer  run-up  distance  for  a  longer 
characteristic  rise  time  at  a  specific  input  pressure.  Therefore,  for  low 
input  pressures  and  critically  long  characteristic  rise  times,  detonation  is 
not  expected  to  occur  in  a  10  cm  bed,  eliminating  the  hazard  of  DDT  in  the 
Case  2  configuration.  In  regards  to  the  actual  rocket  motor  configuration, 
this  implies  that  the  DDT  Case  2  hazard  does  exist  if  the  granular  particles 
are  sufficiently  small  and  the  Zone  1  bed  length  is  greater  than  L$gy  for  the 
corresponding  dP*/dt*. 

Figure  9  shows  the  run-up  to  detonation  distance,  L^gy  as  a  function  of 
t*,  the  compression  wave  rise  time,  for  two  values  of  initial  porosity,  aQ. 
Here,  the  sample  with  the  lower  aQ  is  observed  to  have  a  shorter  run-up  dis¬ 
tance  than  the  more  porous  sample.  This  is  explained  as  follows.  Although 
the  higher  porosity  material  has  a  larger  hot  spot  volume  (as  predicted  by 
Eq.  (16)),  the  low  porosity  material  has  a  much  higher  local  density  and  thus 
contributes  more  energy  to  the  shock  front  upon  decomposition. 
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Figure  8  Prediction  of  run-distance  to  detonation  vs.  input  maximum 
pressure  for  rise-times  increasing  from  1  to  10  usee.  The 
experimental  data  of  J.  J.  Dick  is  shown  as  the  solid  line. 
Both  model  and  data  were  for  a  =  1.5322. 
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CHAPTER  2 

COMPARISON  WITH  EXPERIMENTS 


Introduction 

Clearly  there  is  a  need  to  verify  directly,  or  indirectly,  the  results 
which  the  model  predicts.  One  must  keep  in  mind  that  we  predict  the  run-up  to 
detonation  (time  and  distance)  for  porous  nitramine  explosives  (HMX,  in  par¬ 
ticular)  caused  by  rapid  compression.  As  was  defined  in  the  previous  chapter, 
by  "rapid  compression"  we  mean  rise  rates  to  peak  pressures  of  the  order  of 
1  GPa  within  times  of  many  microseconds. 

The  variables  affecting  the  run-up  conditions  which  are  most  important 
(besides  material  properties)  are: 

(a)  Porosity,  *  ;  or  a  =  l/(l-<t>0);  (typically  1.1  <  a  <  1.5) 

(b)  Peak  pressure,  P*;  (typically  0.1  to  10  GPa) 

(c)  Pressure  rise  time  to  P*,  defined  as  t*  (typically  1  to  10  ysec). 

Review  of  the  open  literature  which  deals  with  the  initiation  to  detona¬ 
tion  of  porous  explosives  indicates  that  although  some  data  exists  for  the 
shock  initiation  (t*  <  1  usee)  of  porous  explosives,  for  our  range  of  a  and  P* 
there  is  no  data  for  ramp  compression  processes,  i.e.  t*  >  1  usee. 

Nevertheless,  we  believe  that  the  model  should  be  tested  against  the  ex¬ 
isting  shock  initiation  experiments.  If  the  predictions  match  the  observa¬ 
tions,  we  will  have  increased  confidence  in  the  usefulness  of  the  model.  Con¬ 
versely,  if  the  data  shows  that  the  model  insufficiently  describes  the  run-up 
conditions,  we  will  be  further  motivated  to  alter  some  of  the  assumptions  and 
the  subsequent  analyses  that  go  into  it. 

Shock  Initiation  Experiments  [Porous  Explosives] 

Dick  [13]  performed  shock  initiation  experiments  using  75  mm  diameter  cy¬ 
lindrical  samples  of  pour  density  (1.24  +  0.04  g/cc,  a  =  1.5322)  HMX.  In  his 
experimental  setup,  the  sample  was  attached  to  an  aluminized  PMMA  driver 
plate,  and  a  rotating-mirror  camera  was  used  to  record  the  shock  transit  time 
through  the  explosive.  Dick  obtained  his  data  points  by  utilizing  a  statis¬ 
tical  analysis  of  large  amounts  of  data  to  estimate  mean  material  response. 


It  was  possible  to  determine  run  distance  and  time  to  detonation  for  a  given 
input  stress  through  transit  time  measurements  at  a  number  of  sample  thick¬ 
nesses. 

Several  calculations  were  made  using  an  initial  porosity  (a  =  1.532) 
which  matched  the  one  used  in  experiments  by  Dick  [131.  The  results  a.  e  shown 
in  Fig.  1G,  a  "Pop-plot"  of  P*  versus  1-sDT*  As  figure  indicates,  the 
slopes  agree  well  for  P*  >  2.2  GPa.  For  1.3  GPa  <  P*  <  2.2  GPa,  our  calcula¬ 
tions  show  little  change  in  run-up  distance  as  a  function  of  P*.  As  P*  drops 
below  1.3  GPa,  1<jqj  is  predicted  to  rise  sharply.  In  fact,  our  code  failed  to 
show  a  detonation  occurring  up  to  t  =  35  usee  if  P*  <  1.1  GPa.  This  pre¬ 
diction  does  not  match  Dick's  data  in  this  region,  indicating  that  our  model 
will  need  to  be  changed  in  order  to  match  the  "low"  pressure  initiation. 

However,  it  goes  without  question  that  the  model  does  a  surprisingly  ex¬ 
cellent  job  in  predicting  the  sensitivity  of  the  run-to-detonation  distance  as 
a  function  of  the  shock  pressure.  No  adjustments  in  any  of  the  property  con¬ 
stants  or  the  kinetic  rates  for  hot-spot  decomposition  were  necessary. 

From  Fig.  10,  one  can  fit  a  relation  for  both  data  and  predictions.  The 
data  of  Ref.  [131  gives  approximately  that  lS0T  *  4.6/P0*56  while  the  model 
gives  LSDT  -  1.71/P0*49  if  P*  >  2  GPA. 

Lindstrom  [14]  examined  the  plane-wave  shock  initiation  of  a  large  part¬ 
icle  size  tetryl  explosive  at  several  pressing  densities,  ranging  from  essen¬ 
tially  TMO  (a  -  1)  down  to  75.1  percent  TMD  (a  =  1.332).  Using  the  explosive 
wedge  technique,  the  shock  buildup  to  detonation  was  observed  and  recorded 
with  a  rotating  mirror  camera.  From  this  film  record,  the  run-distance  and 
time  to  detonation  were  obtained. 

Data  was  obtained  for  several  initial  shock  strengths,  which  were  dic¬ 
tated  by  the  type  of  attenuator  used.  The  inert  attenuator  consisted  of 
either  a  single  metal  driver  plate  or  two  metal  plates  separated  by  a 
liquid.  The  driver  plate  free  surface  velocity  was  computed  from  charged-pin 
detector  data  and  from  this  other  shock  parameters,  such  as  Initial  shock 
strength,  shock  velocity,  and  particle  velocity,  were  derived. 

Figure  11  Is  the  "Pop-plot"  comparing  results  from  our  calculations  to 
experimental  results  obtained  by  Lindstrom  [14J  for  a  porosity  of  s  *  1.332. 
Note  that  our  calculations  are  representative  of  HMX  whereas  Lindstrom' s  data 
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Figure  10  A  closer  view  of  the  predictions  shown  In  Figure  8  (run-distance 
vs.  P*).  Comparison  of  the  data  by  J.  J.  Dick  with  predictions 
assuming  t*  =  1  usee  show  good  agreement. 
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distance  to  detonation  vs.  maximum  input  pressure  for 
1.332.  The  solid  line  represents  data  by  I.  E.  Lindstrom 
tetryl.  The  dash  line  represents  predictions  for  HMX  with 
same  porosity  using  a  pressure  rise  time  of  t*  =  1  nsec. 


is  for  tetryl.  Nonetheless,  the  two  are  in  very  close  agreement  for 
P*  >  2  GPa.  As  P*  drops  below  2  GPa,  however,  our  calculated  run-up  distances 
exhibit  unexpected  behavior.  There  is  a  minimum  at  P*  =  1.5  GPa  and  below 
this  value  our  L^y's  rise  sharply.  The  positive  slope  for  1.5  GPa  <  P*  <  2 
GPa  may  not  be  physically  meaningful  since  increasing  the  shock  strength 
should  not  increase  the  run-up  distance.  Again,  we  will  need  to  consider  what 
aspects  of  the  shock  initiation  process  at  the  relatively  low  pressures  are 
not  being  described  by  our  model. 

Sandusky  and  Bernecker  [151  obtained  shock  initiation  data  for  several 
porous  explosives,  including  class  0  HMX  with  a  porosity  of  a  =  1.37  and  class 
E  HMX  with  a  porosity  of  a  =  2.27.  In  their  experiments,  the  explosive  was 
radially  confined  and  subjected  to  a  long  duration  (>20C  us),  low-level  shock 
loading  (<0.18  GPa). 

The  long  duration  loading  was  accomplished  by  using  a  powder  gun  to  pro¬ 
pel  a  25.4  mm  diameter  by  305  mm  long  piston  into  the  tube  that  contained  the 
porous  bed.  It  was  generally  found  that  when  high  confinement  steel  tubes 
were  used,  the  reaction  grew  rapidly  and  built  to  detonation.  Often  for  these 
experiments,  the  inner  wall  of  the  tube  fragments  showed  the  location  of  the 
onset  of  detonation.  In  low  confinement  tubes,  higher  loadings  were  required 
for  a  transition  to  detonation. 

Typical  run-up  distances  were  large  at  these  low  initiation  shock  pres¬ 
sures.  For  example,  for  HMX  with  a  =  1.37,  Sandusky  and  Bernecker  reported 
[ 15]  the  fol lowing: 

LSDT 
39  mm 

50  mm 

An  alternative  to  the  run  distance  versus  pressure  plot,  the  so-called 
"Pop-plot,"  in  correlating  shock  initiation  data  is  a  log-log  plot  of  run-up 
distance  versus  the  particle  (material)  velocity.  Sandusky  and  Bernecker  pre¬ 
sented  a  composite  plot  of  L<jgy  versus  up  containing  data  for  several  differ¬ 
ent  explosives  for  a  variety  of  porosities,  including  nearly  voidless  explo¬ 
sives.  They  showed  that  almost  all  the  data  fit  within  a  band  of  constant 
slope.  A  straight-line  fit  to  the  data  gives  the  relation  L^qj  »  1.6  Up"^ 


Shock  Pressure,  P* 
0.108  GPa 
0.1732  GPa 
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where  L$0j  is  in  mm  and  up  is  in  mm/usec.  The  data  covers  particle  velocities 
ranging  from  approximate^  0.1  to  1.2  mm/usec.  This  regime  includes  the  very 
low  shock  initiation  pressures  such  as  those  used  by  Sandusky.  We  believe 
this  type  of  correlation  is  a  slight  improvement  tc  that  of  the  Pop-plot  since 
a  wider  spectrum  of  data  falls  within  a  constant  slope  region.  In  particular, 
it  permits  the  inclusion  of  data  using  low  shock  initiation  pressures  and  data 
from  nearly  homogeneous  explosives. 

In  Figs.  12  and  13,  results  from  our  calculations  for  a  =  1.5322  are 
shown  as  a  function  of  particle  velocity  behind  the  shock  front.  Figure  12 
shows  run-up  distance  to  detonation,  and  Fig.  13  shows  time  to  detonation. 
Here,  up  was  calculated  at  the  left  boundary  of  the  cast  explosive,  adjacent 
to  the  granulated  bed.  The  solid  line  in  each  figure  is  an  "eyeball"  best  fit 
of  our  calculations,  which  are  shown  by  the  solid  circles.  The  results  are 
qualitatively  correct  in  that  both  run-up  distance  and  time  to  detonation  de¬ 
crease  linearly  with  increasing  shock  strength  (represented  by  up).  More 
quantitative  data  for  this  porosity  (a  =  1.5322)  is  needed  in  the  future  to 
enhance  this  parametric  study. 

Figures  14  and  15  repeat  the  Ljqj  versus  up  plots  for  an  Initial  porosity 
of  a  =«  1.332.  This  porosity  matches  the  one  used  by  Lindstrom,  whose  data  for 
porous  tetryl  is  also  included  in  these  two  figures.  In  each  figure,  the  open 
circles  are  lindstrom's  data  points  and  the  solid  line  is  an  "eyeball"  best 
fit  of  this  data.  The  solid  circles  are  results  from  our  computer  code.  In 
Fig.  14,  tnree  of  our  calculations  correlate  well  with  the  data.  Two  predic¬ 
tion  points,  well  to  the  left  of  the  data  line,  are  results  from  runs  using 
P*'s  of  1.2  and  1.5  GPa.  These  points  are  from  the  troublesome  region  of 
Fig.  11,  discussed  earlier.  The  dashed  line  in  Fig.  15  is  an  "eyeball"  best 
fit  of  our  calculations.  Even  though  the  slope  was  intended  to  match  that  of 
Lindstrom's  data,  our  calculations  still  show  good  qualitative  agreement. 
Notice  that  there  are  no  "stray"  points  in  the  time  to  detonation  plot  as 
there  were  in  the  run-up  distance  to  detonation  plot  in  Fig.  14.  In  other 
words,  as  the  initial  peak  input  pressure  was  varied,  the  time  to  detonation 
did  not  exhibit  the  unexpected  behavior  that  was  found  in  the  run-up  distance 
to  detonation. 

All  of  the  data  obtained  by  Sandusky  for  porous  HMX  (as  well  as  for  other 
explosives)  was  for  initiation  by  shocks  of  strengths  less  than  0.18  GPa, 


T 


T 


Detonation  Run  Distance  vs  Particle  Velocity 

a  =  1.5322  (Static  Pore  Collapse} 
t*  =  1.0/xsec  - Best  Fit 


p*  =  1.2  GPa 


p*  =  1.5  GPa 


p*  =2.0  GPa 


p4=3.0GPa- 


0.2  0.3  04  05  0.6  0.7  0.80.910 

Particle  Velocity  (mm//xsec) 


Figure  12  Predicted  run-distance  vs.  imposed  particle  velocity  at  the 
Zone  1/Zone  2  interface.  (The  peak  input  pressures  needed  to 
provide  the  particle  velocities  are  indicated  in  the  figure.) 


Data  (Ref  14) 


Detonation  Run  Distance  vs 
Particle  Velocity 

a  =  1.332  Static  Pore  Collapse 
t*=  l.Oyu-sec  •  DSDT  Code 

0  Lindstrom  ( Tetryl ) 


0.2  0,3  ,0.4  0.5  0.6  0.7  0.80.910 

Particle  Velocity  (mm/yusec) 


Figure  14  Predicted  run-distance  vs.  particle  velocity 
compared  with  the  data  of  I.  E.  Lindstrom. 


which  resulted  in  very  large  run-up  distances.  In  order  for  us  to  compare  our 
calculations  to  this  data,  it  is  necessary  to  expand  the  time  limits  set  on 
our  computer  runs.  Even  though  this  will  likely  be  very  expensive,  limited 
runs  of  this  nature  will  be  done  in  the  future.  As  previously  pointed  out, 
our  model  does  not  do  as  good  a  job  in  matching  experiments  as  P*,  the  peak 
input  pressure,  drops  below  2  GPa.  By  expanding  the  computer  time  and  carry¬ 
ing  out  more  low  input  pressure  calculations,  it  is  hoped  that  trends  will  be 
found  that  provide  a  clue  to  the  shortcomings  of  the  model  in  this  region. 
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CHAPTER  3 

MODEL  EVALUATION:  SENSITIVITY  STUDIES 
ON  PORE  COLLAPSE  (ANO  OTHER  COMPONENTS) 

Arrhenius  Kinetics 

The  first  sensitivity  study  presented  here  is  that  of  the  "activation 
temperature"  term  which  Is  found  in  Eq.  (1.13).  "Activation  temperature"  is 

the  name  given  to  the  variable  grouping  E*/R,  where  E*  is  the  activation 

energy  and  R  Is  the  universal  gas  constant.  Physically,  as  used  in  the 
Arrhenius  burn  model  of  Eq.  (1.13),  the  activation  energy  represents  the 

minimum  value  of  molecular  energy  needed  to  guarantee  chemical  reaction. 
Figure  16  is  a  plot  of  run-up  distance  to  detonation  as  a  function  of 

activation  temperature,  or  ACTEMP,  for  two  different  values  of  initial 
porosity.  Here,  ACTEMP  is  normalized  with  respect  to  the  "standard"  ACTEMP  of 
14400  K,  which  is  the  value  used  for  all  other  calculations  presented  in  this 
report.  A  ramp  wave  input  with  a  peak  pressure  of  2  GPa  and  a  rise  time  of 
10  nsec  was  used  for  these  runs,  as  well  as  the  static  pore  collapse  model  (to 
be  discussed  in  the  next  section)  and  a  finite  difference  grid  spacing 
of  ax  *  0.5  mm. 

Upon  examining  the  figure,  both  curves  initially  show  steep  increases  in 
L^qj  as  ACTEMP  increases.  Intuitively,  this  is  expected  since  raising  the  ac¬ 
tivation  energy  means  more  energy  is  required  to  initiate  reaction.  More  time 
is  needed  for  the  coalescing  stress  waves  to  build  up  to  the  required  energy 
level.  During  this  time  the  stress  front  continues  to  propagate  through  the 
explosive.  Thus,  reaction  is  finally  initiated  further  downstream.  Surpris¬ 
ingly,  the  similarities  between  the  two  curves  end  after  their  initial  as¬ 
cents.  The  solid  curve,  for  a  =  1.1176,  levels  out  only  temporarily  before 
sharply  rising  again,  showing  no  further  indication  of  reaching  a  steady 
value.  On  the  other  hand,  the  dashed  curve,  for  a  *  1.2667,  reaches  a  peak 
before  dropping  slightly  to  a  constant  value.  Further  Investigation  is  needed 
to  clarify  why  the  more  porous  sample  is  not  affected  by  changes  in  ACTEMP  in 
this  constant  region.  Clearly  other  mechanisms  are  involved. 


Predicted  run-distance  to  detonation  as  a  fu 
of  the  activation  energy  for  explosive  decom 
Two  different  initial  porosities  are  studied 
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Pore  Collapse  Model:  Static  Compaction  Assumption 

As  mentioned  in  Chapter  1,  the  pore  collapse  process  is  represented  by  a 
three  step  static  compaction  model.  A  detailed  derivation  of  this  model  is 
given  by  Carroll  and  Holt  in  Ref.  1 8 1 .  A  brief  description  is  presented 
here.  The  static  model  involves  no  rate  dependence  term  in  describing  the 
pore  collapse  process.  The  relation  used  is 

P  •  Peq(.)  (3.1) 

where  P  is  the  applied  pressure  and  P0q  is  the  equilibrium  pressure  of  the 
solid  corresponding  to  the  instantaneous  porosity.  For  instantaneous  material 
response  (an  assumption  inherent  with  the  use  of  the  static  model),  P  is  also 
equal  to  the  pressure  of  the  aggregate  porous  material,  which  is  a  volume 
average  of  the  pressure  of  the  solid  material  calculated  from  an  appropriate 
constitutive  relation  and  the  pressure  of  the  gas  inside  the  pore  (assumed  to 
be  at  atmospheric  pressure). 

Carroll  and  Holt  ( 8 1  define  three  regimes  of  compression:  the  elastic, 
elastic-plastic,  and  the  plastic  phases.  The  equilibrium  pressure  for  each 
phase  and  the  appropriate  range  over  which  each  applies  are  given  by 


elastic  phase  aQ  >  a  a 
4G( a  -a) 

Peq(a)  =  3a(a- 1 ) 


elastic-plastic  phase  >  a  > 

Peq(a)  “  f  Y  t1  “  H  (ao'a)  +  ln 
plastic  phase  a2  >  a  >  1 


2G(a  -a) 

1-y^TT 


1} 


(a)  =  |  Y  t n(r2-r) 


where  the  limits  between  the  phases  are  expressed  as 

2G  aQ  ♦  Y  2G  co 

“1  *  2G  +  Y  and  “2  *  2G  +  Y* 


(3.2) 


(3.3) 


(3.4) 


(3. 5, 3. 6) 


27 


To  calculate  the  porosity  at  each  timestep  in  the  computer  program,  an  iter¬ 
ative  process  is  used  to  solve  Eg.  (1)  utilizing  the  appropriate  relation  for 
Peq.  The  individual  terms  in  the  above  equations  are  defined  in  Table  2. 

Pore  Collapse  Model:  Dynamic  Compaction  Assumption 

Recently  the  static  pore  collapse  model  was  replaced  by  two  different 
versions  of  dynamic  material  response  models.  The  first  version  discussed 
here  was  derived  by  Carroll  and  Holt  (81,  whose  corresponding  static  model  was 
used  previously  in  our  code.  In  fact,  their  static  model  <s  a  direct  result 
of  their  dynamic  model  upon  the  assumption  that  the  applied  pressure  rate  is 
infinitesimally  small.  In  this  section,  Carroll  and  Holt's  dynamic  model  will 
be  investigated  extensively,  including  many  comparisons  and  references  to  the 
static  model.  Following  this,  a  brief  look  at  a  second  dynamic  model,  that  of 
Kooker  and  Anderson  [161,  will  be  made. 


Oynamic  Model  Version  1:  Carroll  and  Holt  [81 

A  detailed  derivation  of  this  model  is  given  in  Ref.  [8 J .  A  brief  des¬ 
cription  is  presented  here.  The  dynamic  relation  is  given  by 

T2  Y  Q(a ,  a,  a)  =  P  -  Pgq(a).  (3.7) 

The  strong  rate  dependence  of  the  pore  collapse  process  is  embedded  in  the 
term  Q(a,  a,  a),  which  is  given  by 

Q(a,  a,  a)  =  -  a( (a-l)'l/3  -  a'1/3l  +  |  a2((a-I)'4/3  -  a’473]  (3.8) 


The  term  t  in  Eq.  (3.7)  is  a  material  time  constant 
dimension  of  time)  defined  by 
a2 

2  *  °  ao 

3Y(ao-I)2/3‘ 

The  terms  in  this  equation  are  defined  in  Table  2. 


(having  the  physical 


(3.9) 


P  and  Peq  are  as  defined  for  the  static  model.  Note  that  P  is  actually 
the  applied  pressure  which,  because  of  the  finite  material  response  times  as- 
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Table  2  Nomenclature 


Symbol 

Oef Initlon 

Units 

Typical  Value 

ao 

initial  porosity 

1.1176 

a 

Instantaneous  porosity 

ao  ~  *  1 

Y 

yield  strength 

MPa 

51.7  (HMX) 

G 

shear  modulus 

MPa 

3516  (HMX) 

ADDITIONAL  TERMS  USED  IN 

DYNAMIC  MODEL 

p 

Initial  density  of  the  solid 

g/cc 

1.90  (HMX) 

ao 

Initial  pore  radius 

•*m 

20 

A 

time  Increment 

nsec 

5 

/ 


"T.  ’  *%  *W 
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sociated  with  the  dynamic  model,  is  not  equal  to  the  pressure  of  the  aggregate 
material,  except  possibly  for  very  small  applied  pressure  rates.  Since  the 
applied  pressure  is  only  known  at  one  location  (and  thus  for  only  one  finite 
difference  cell),  the  material  boundary,  we  are  forced  to  treat  this  term  just 
as  was  done  for  the  static  model.  That  is,  P  is  assumed  to  be  equal  to  the 
pressure  of  the  aggregate  material. 

Upon  examining  Eq.  (3.6),  it  is  evident  that  if  the  ratio  of  yield 
strength  to  shear  modulus  is  small,  which  is  generally  the  case  (for  HMX, 
Y/G  *  0.0147),  the  change  in  porosity  during  the  elastic  and  elastic-plastic 
K.iases  is  small.  Thus,  to  simplify  the  analysis  of  Eq.  (3.7)  a  good 
approximation  may  be  made  by  neglecting  rate  dependence  during  these  first  two 
phases.  In  doing  so,  Eq.  (3.7)  can  be  rewritten  by  substituting  Eq.  (3.4) 
for  Peq(a), 

t2  Y  Q(a,  i,  a)  -  P  -  |  Y  in  (yyy).  (3.10) 

Carroll  and  Holt  [8|  used  the  following  finite  difference  form  of 
Eq.  (3.10)  for  numerical  analysis. 

aN*l  3  2aN  *  Vl  +  11/fl(aM),(i  (aN-aN-l)Z  f4(aN} 

-  (a2/'2Y)(Pn  -  \  Y  in  (3.11) 

H 

Where  f$(a)  =  (a-i)'5/3  -  a‘S/3. 

In  summary,  the  dynamic  model  consists  first  of  modeling  statically  the 
elastic  and  elastic-plastic  phases  of  compaction,  whpre  the  change  in  porosity 
is  very  small.  When  the  transition  porosity  a.,  is  attained,  Eq.  (3.11)  is 
used  to  calculate  new  values  of  porosity.  The  ini  dal  conditions  needed  for 
Eq.  (3.11),  which  are  and  are  those  values  calculated  from  the  pre¬ 

vious  two  timesteps.  Note  that  this  form  of  the  dynamic  model  gives  an  ex¬ 
plicit  expression  for  a.  No  Iteration  is  required. 

The  first  comparison  between  the  static  and  dynamic  models  will  be  that 
of  run-up  distance  to  detonation  ( L^qj)  as  a  function  of  initial  porosity. 
The  results  are  plotted  in  Fig.  17.  For  these  runs,  an  input  ramp  wave  with  a 
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peak  pressure  of  2  GPa  and  a  rise  time  of  10  usee  was  used,  along  with  a 
finite  difference  grid  spacing  of  0.5  mm.  Runs  using  the  dynamic  model  also 
required  the  initial  pore  radius  (aQ)  to  be  Input.  Here,  20  nm  was  used. 

Upon  examining  Fig.  17,  the  two  models  predict  the  same  general  trend 
over  the  whole  range  of  porosities.  The  curves  initially  show  rapid  declines 
in  L$Dy  since  the  more  porous  materials  have  a  higher  hot  spot  volume,  and 
thus  more  irreversible  energy  is  deposited  to  initiate  reaction.  If  the  two 
curves  were  extended  to  include  lower  porosities,  each  would  approach  an  asym¬ 
ptote  at  around  a  *  1.02.  This  means  that  for  porosities  at  or  below  this 
value,  detonation  ceases  to  occur  in  our  10  cm  bed  for  the  given  ramp  input. 
However,  for  a  >  1.35,  the  Logy's  start  to  rise.  In  this  range  of  higher  por¬ 
osities,  even  though  the  hot  spot  volumes  are  greater,  the  local  density  of 
the  material  has  decreased  to  a  po^nt  where  its  effect  overrides  that  of  hot 
spot  volume.  As  previously  concluded  for  the  static  pore  collapse  runs,  the 
higher  porosity  material  (i.e.  the  material  with  the  lower  local  density)  has 
less  mass  and  thus  less  energy  upon  decomposition  available  to  support  a  shock 
front.  Thus  the  distance  at  which  a  steady  state  detonation  is  reached  Is  in¬ 
creased.  Although  just  a  conjecture  here,  the  higher  porosity  material  is 
more  likely  to  bend  and  distort  the  shock.  This  possibility  could  be  studied 
if  the  present  analysis  were  extended  to  two  dimensions. 

Another  factor  to  discuss  is  the  equivalent  volume  concept  inherent  in 
the  hot  spot  model  used  here;  that  is,  the  initial  volume  of  hot  spot  material 
is  equal  to  the  original  pore  volume.  This  assumption  is  reasonable  for,  and 
intended  for,  low  porosities.  This  fact  should  be  kept  in  mind  when  consider¬ 
ing  higher  porosities  where  the  concept  may  not  be  accurate.  A  revision  of 
the  equivalent  volume  concept  for  high  porosities  (a  >  1.25),  which  is  a  topic 
of  current  investigation,  might  significantly  affect  the  results. 

Referring  back  to  Fig.  17,  for  a  >  1.19,  the  dynamic  model  predicts 
lSDTs  that  are  at  or  be1ow  those  predicted  by  the  static  model  while  the  re¬ 
verse  Is  true  for  a  <  1.19.  In  addition,  the  Logy's  from  the  dynamic  model 
are  less  sensitive  to  changes  In  porosity  in  the  range  1.25  <  a  <  1.53.  The 
static  curve  has  a  well  defined  minimum  at  a  *  1.36  whereas  the  dynamic  curve 
has  a  less  apparent  minimum  at  a  *  1.33. 

Table  3  compares  Ljgy  and  t$gy  (time  to  detonation)  values  from  the 
static  and  dynamic  models  for  a  shock  wave  Input  (t*  *  1  usee)  as  well  as  for 


=2.0  GPa 


Figure  17  Predicted  run-distance  to  detonation  vs.  porosity 

comparing  the  static  pore  collapse  assumption  to  the 
dynamic  pore  collapse  model  of  Carroll  and  Holt. 


Table  3  Comparison  of  Run  Up  Distance  and  Time  to  Oetonatlon 
a  -  1.2667,  P*  *  2,  t  *  1,  10 


ramp  p  (mm) 

shock 

t*  =  10  usee 

t*  =  1  usee 

Static 

13.4728 

1.972 

Oynamic 
(aQ  »  20  um) 

12.9756 

1.9495 

X  Olff 
(of  dynamic) 

-3.7X 

-1.1X 

t soT  (•iSec) 

t*  =  10  usee 

t*  =  1  usee 

Static 

11.23 

2.095 

Oynamic 
(aQ  =»  20  um) 

10.99 

2.05 

X  Olff 
(of  dynamic) 

-2. IX 

-2.  IX 
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the  ramp  wave  input  just  discussed.  The  porosity  considered  here 
is  a  =  1.2667.  This  comparison  is  done  to  show  that  the  dynamic  compaction 
model  gives  consistent  results  (compared  to  the  static  model)  for  input  rise 
times  approaching  a  true  shock.  The  consistency  of  the  static  model  for  these 
inputs  has  already  been  demonstrated. 

The  variation  of  L^gy  with  peak  pressure  input  for  both  compaction  models 
is  shown  in  the  "Pop  plot"  in  Fig.  18  for  a  =  1.5322.  The  solid  line  on  the 
plot  is  experimental  data  obtained  by  Dick  [ 13 j  for  this  same  porosity.  The 
results  clearly  Indicate  almost  no  deviation  between  the  static  and  dynamic 
models  for  the  shock  input  (t*  =  1  usee)  and  only  a  slight  deviation  when  P*  < 
2.5  GPa  for  the  ramp  input  (t*  =  10  usee). 

The  next  parametric  study  examines  the  effect  of  initial  pore  radius  on 
*-SDT  fr°r  the  dynamic  compaction  model.  Figure  19  is  a  plot  of  L$gy  versus  aQ 
for  four  different  porosities  using  a  ramp  input  of  P*  =  2  GPa  and 
t*  =  10  usee.  For  aQ  <  40  urn,  the  run  up  length  in  each  porous  sample  is 
virtually  independent  of  aQ.  However  at  aQ  =  40  pm,  L^gy  in  the  a  =  1.2667 
sample  drastically  drops.  This  is  assumed  to  be  a  "trouble  spot"  with  this 
dynamic  model.  For  aQ  >  40  um,  the  a  =  1.2667  and  a  =  1.1875  curves  oscillate 
as  aQ  increases,  which  also  can  be  attributed  to  numerical  difficulties. 

Further  evidence  of  numerical  difficulties  with  the  dynamic  model  is 
shown  in  Fig.  20,  which  shows  the  hot  spot  temperature  profile  about  2.5  usee 
after  detonation  for  a  grid  spacing  of  0.5  mm,  a  porosity  of  1.1176,  a  pore 
radius  of  39  um,  and  a  ramp  input  of  P*  =  2  GPa  and  t*  =  10  usee.  The  profile 
is  marked  by  sharp  oscillations,  which  are  a  result  of  the  explicit  finite 
difference  relationship  used  to  calculate  the  porosity  at  each  timestep. 

Figure  21  is  a  hot  spot  temperature  profile  using  the  static  compaction 
model  with  all  inputs  being  the  same  as  those  for  Fig.  20.  Here  the  oscilla¬ 
tions  are  very  small,  almost  negligible,  indicating  a  greater  stability  of  the 
static  model.  Note  also  that  the  average  hot  spot  temperature  is  about  800  K 
greater  using  the  dynamic  model. 

Next,  the  effect  of  grid  size  and  the  length  of  the  timestep  on  lSDT  are 
examined.  The  finite  difference  cell  size  used  by  Krier  and  Cudak  [lj  and 
previously  in  this  research  was  0.5  mm.  However,  reducing  the  cell  size  to 
0.25  mm  caused  Ljqj  to  go  up  by  about  10  percent  regardless  of  whether  the 
static  or  dynamic  model  was  used.  Run  up  lengths  from  the  static  model  were 


35.0 


Figure  18  Run-distance  to  detonation  vs.  peak  input  pressure 

comparing  data  by  J.  J.  Dick  with  predictions  for  both 
dynamic  and  static  pore  collapse. 
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Figure  21  Hot  spot  temperature  profile  in  the  explosive  bed  about  2.5 
after  detonation  using  the  static  po-e  collapse  model. 
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virtually  Insensitive  to  further  reductions  In  cell  size  down 
to  ax  =  0.125  mm.  On  the  other  hand,  run  up  lengths  from  the  dynamic  model 
dropped  by  more  than  63  percent  when  0.167  mm  <  ax  <  0.20  mm!  These  results 
are  summarized  In  Table  4. 

Figures  22  and  23  present  hot  spot  temperature  profiles  for  runs  Ident¬ 
ical  to  those  from  Figs.  20  and  21,  respectively,  except  that  ax  has  been  re¬ 
duced  from  0.5  mm  to  0.25  mm.  (The  profiles  are  narrower  in  Figs.  22  and  23 
because  detonation  occurred  later  for  ax  =  0.25  mm.)  Comparing  the  two  static 
runs,  the  oscillations  in  Fig.  23  are  sharper,  yet  they  remain  very  small  and 
stable.  Also,  the  values  of  the  hot  spot  temperatures  for  the  two  runs  are 
only  slightly  different.  Comparing  the  twe  dynamic  runs,  the  oscillations  in 
Fig.  22  remain  large,  as  they  were  in  Fig.  20.  Furthermore,  the  reduction 
in  ax  caused  an  increase  in  the  average  hot  spot  temperature  by  more  than 
2000  K! 

The  timestep  used  in  all  previous  runs  using  the  static  model  was 
5  nsec.  Extending  this  analysis  to  Include  timesteps  ranging  from  4  nsec  to 
10  nsec  gives  run  up  lengths  that  do  not  vary  significantly.  However,  the 
dynamic  model  gives  LSDj * s  varying  by  up  to  11  percent  in  the  range  of 
4  nsec  <  at  <  10  nsec.  These  results  are  summarized  in  Table  5. 

It  was  initially  conjectured  that  the  numerical  difficulties  with  the 
dynamic  model  could  be  due  to  an  incorrect  artificial  viscosity,  a  term  which 
is  included  in  the  momentum  equation  to  account  for  irreversibilities  across 
the  shock  wave.  However,  the  sharp  oscillations  in  the  hot  spot  temperature 
remained  for  all  values  of  artificial  viscosity  tried.  Figure  24  is  a  plot  of 
L$0t  versus  ARV,  the  coefficient  of  quadratic  artificial  viscosity.  The  value 
of  ARV  used  in  all  runs  up  until  now  was  3.0.  Included  here  are  the  static 
results  as  well  as  dynamic  results  for  pore  radii  of  39  urn  and  50  urn.  Other 
inputs  were  a  =  1.1176,  P*  =  2  GPa,  t*  =  10  usee,  and  ax  =  0.25  mm.  The 
static  results  given  by  the  solid  line  show  a  small,  smooth  increase  in  Lso-j- 
as  ARV  increases.  The  dynamic  results,  however,  once  again  exhibit  very  er¬ 
ratic  behavior.  There  seem  to  be  no  similarities  between  the  two  dynamic 
curves  representing  the  two  pore  radii.  Recall  that  Fig.  19,  where  ax  = 
0,50  mm,  showed  no  dependence  of  l$Qy  on  pore  radius  for  a  =  1.1176, 
P*  =  2  GPa,  t*  =  10  nsec.  Note  also  the  "trouble  spot"  at  ARV  =  2.0  for 
a0  =  39  »m.  To  further  announce  the  instability  of  the  dynamic  method,  this 
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Table  4  Effect  of  Finite  Difference  Cell  Size  on  Computed  Results  Using 
(a)  The  Dynamic  Model  and  (b)  The  Static  Model  (P*  =  2  GPa,  t* 

=  10  usee,  a  =  1.1176,  aQ  =  39  urn,  At  3  5  nsec) 


Run 

ax  (mm) 

L$dt  (mni) 

Deviation  of 

L$qt  from  Run  1 

tSDT  (ySec^ 

1 

0.50 

20.9486 

12.45 

2 

0.25 

23.0206 

+  10% 

13.22 

3 

0.20 

7.5382 

-  64% 

8.17 

4 

0.167 

7.6951 

-  63% 

(a) 

6.26 

Run 

ax  (mm) 

lsqt  (mm) 

Deviation  of 

L$dt  from  Run  1 

t$QT  (usee) 

1 

0.50 

20.5049 

12.41 

2 

0.25 

22.5351 

+  lose 

13.13 

3 

0.20 

22.8153 

+  11% 

13.155 

4 

0.167 

23.0246 

+  12% 

13.235 

5 

0.125 

23.1690 

+  13% 

13.280 

(b) 


=  1.1176 
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Table  5 

Effect  of  Timestep  on  Computed  Results  Using  (a) 

The  Dynamic  Model 

and  (b)  The  Static  Model  (P*  * 

2  GPa,  t*  «  10  usee,  a  -  1.1176.  a- 

*  39  urn,  ax 

*  0.25  mm) 

U 

Run 

At  (nsec) 

tS0T  (nwi) 

Deviation  of 
lSDT  from  Run  * 

tSDT  (usec) 

1 

5 

23.0206 

13.22 

2 

4 

20.5354 

-  11* 

13.224 

3 

6 

22.5460 

-  2% 

13.158 

4 

10 

21.2851 

-  8X 

12.54 

(b) 

Run 

At  (nsec) 

lS0T  (m) 

Deviation  of 
L$ot  from  ^un  1 

^SDT  (uSec) 

1 

5 

22.5351 

13.13 

2 

4 

22.5272 

-  0.03X 

13.116 

3 

6 

22.5415 

+  Q.03X 

13.14 

4 

10 

22.7974 

+  1.16* 

13.26 

(b) 
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"trouble  spot"  can  be  avoided  simply  by  changing  either  ax,  at,  ARV,  or  aQ 
only  slightly! 

Finally,  Figs.  25  and  26  give,  respectively,  dynamic  and  static  hot  spot 
temperature  profiles  for  a  *  1.5322.  Other  Input  data,  P*  *  2  GPa, 
t*  *  10  usee,  and  ax  *  0.5  mm,  is  Identical  to  that  of  Figs.  20  and  21  except 
that  here  a  pore  radius  of  20  um  Is  used.  Most  noticeable  Is  much  lower 
values  of  the  hot  spot  temperature  for  the  more  porous  material.  In  addition, 
the  oscillations  In  the  dynamic  curve  are  much  smaller,  fewer  In  number,  and 
not  as  sharp  as  they  were  In  Fig.  20  for  the  less  porous  material.  The 
average  value  of  the  hot  spot  temperature  appears  to  be  about  equal  upon  com¬ 
parison  of  the  static  and  dynamic  results,  whereas  Figs.  20  and 
21  (a  *  1.1176)  show  an  800  K  difference  in  the  average  hot  spot  temperature 
between  the  two  models. 

In  summary,  a  dynamic  pore  collapse  model  should  be  a  more  realistic 
representation  of  the  actual  compaction  process.  In  the  presentation  here, 
however,  certain  combinations  of  ARV,  ax,  and  a0  produce  inexplalnable 
"trouble  spots"  in  which  the  run  up  distance  and  hot  spot  temperature  profl’e 
deviate  greatly  from  the  expected  trends.  This  Is  assumed  to  be  a  major  fault 
in  this  dynamic  model.  These  "trouble  spots"  are  not  found  when  using  the 
static  pore  collapse  model. 


Dynamic  Model  Version  2:  Kooker  and  Anderson  (16| 

This  section  briefly  describes  the  dynamic  pore  collapse  model  of  Kooker 
and  Anderson  [  16 ]  as  well  as  how  it  was  incorporated  into  our  DSDT  code.  A 
comparison  is  then  made  of  results  from  this  model  and  Carroll  and  Holt's  dyn¬ 
amic  model. 

The  dynamic  pore  collapse  relation  of  Kooker  and  Anderson  [161  was  writ¬ 
ten  specifically  for  the  modeling  of  hivelite  solid  propellant  combustion. 
Their  theory  is  similar  to  Carroll  and  Holt's  [8 1  in  that  it  is  based  on  Iso¬ 
lated  pores  or  voids.  They  give  the  following  time-dependent  relation  for  a. 


na  D 
a(a-l)  Dt 


(a) 


■psm(oc’a’ec)  +  Peq(a’ec)  +  Pgas(a)* 


(3.12) 


This  Is  the  radial  momentum  equation  through  a  hollow  sphere,  neglecting  In¬ 
ertia.  Psn|  is  the  spherical  stress  in  the  solid  matrix  material,  Pgas  Is  the 


=  1.5322 
=2.0  GPa 
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pressure  of  the  gas  trapped  in  the  pores,  and  Peq  is  again  the  equilibrium 
stress  field  In  the  solid  surrounding  the  pore.  Also  in  Eq.  (3.12),  na  is  a 
viscosity  type  coefficient,  p  is  the  density  of  the  condensed  phase,  and  er 

v*  L 

is  the  internal  energy  of  the  condensed  phase. 

Different  from  Carroll  and  Holt's  model,  Kooker  and  Anderson's  model  does 
not  Involve  separate  phases  of  compression.  They  point  out  that  while  Carroll 
and  Holt's  model  provides  an  excellent  description  for  porous  metals,  it  works 
poorly  for  solid  propellant.  Therefore,  Kooker  and  Anderson  present  a  differ¬ 
ent  expression  for  peq(a*ec)»  which  is 

Peq<’'ec>  '  Vf<“-ec>  tn  <J>  (3-13) 


where 

Teff(a>ec)  =  (ro  +  (TrT0)(1  -  c  1  -  ^2(l-c-1)  1 }  •  (1  -  A-^)2}  (3.14) 


c  h 


a  -  1 


Teff  plays  the  role  of  a  porosity-dependent  yield  stress.  The  multiplier  in¬ 
volving  internal  energy  in  Eq.  (3.14)  provides  for  thermal  softening,  i.e.  at 
the  melting  condition  the  material  is  unable  to  resist  pore  collapse.  eQ  is 
the  initial  value  of  ec,  and  eml  is  the  melt  energy  of  the  condensed  phase- 
above  e0.  The  parameters  (vV^Pj)  in  EcJ*  (3-14)  are  determined  exper¬ 
imental  ly. 


The  following  finite  difference  form  of  Eq.  (3.12)  was  incorporated  into 
Krier  and  CudaK's  (1]  code  for  the  calculation  of  a  at  each  timestep. 


N+l 

a 


p"  ♦  p"  i  ,  ♦  .n 
eq  gas1 


(3.15) 


This  relation  Is  used  throughout  the  entire  compaction  process.  Note  that 
since  inertia  is  neglected,  this  equation  is  first-order  (compare  with  that  of 
Carroll  and  Holt,  Eq.  (3.11),  which  is  second-order).  The  size  of  the  time- 
step,  A,  is  5  nsec. 

The  use  of  this  model  requires  values  for  the  several  parameters 
mentioned  above.  Since  these  values  are  not  available  for  our  explosive,  HMX, 
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we  are  forced  to  use  those  values  given  by  Kooker  and  Anderson  for  hivellte 
1C86-8A,  which  are  listed  in  Table  6.  Even  though  the  use  of  these  values  ob¬ 
viously  will  not  give  accurate  numerical  results  for  HMX,  It  is  the  intention 
here  to  at  least  incorporate  the  model  in  a  working  form  into  the  code.  Ac¬ 
curate  values  for  those  parameters  listed  in  Table  6  are  hoped  to  be  obtained 
in  the  future  for  HMX. 

Finally,  keeping  in  mind  the  inaccuracies  mentioned  above,  the  results 
from  this  model  will  be  compared  to  results  from  the  Carroll  and  Holt  model 
for  a  ramp  wave  input  of  P*  ■  2  GPa,  t*  3  10  usee,  an  initial  porosity 
of  aQ  *  1.1176,  and  a  finite  difference  grid  spacing  of  ax  =  0.5  mm.  The  com¬ 
parison  is  summarized  in  Table  7.  Examining  this  table,  it  is  seen  that  while 
the  run  up  distance  given  by  Kooker  and  Anderson's  model  is  25  percent  lower 
than  the  Ljqj  from  Carroll  and  Holt's  model,  the  time  to  detonation  given  by 
the  former  is  extremely  high  --  142  percent  higher  than  the  t^y  from  Carroll 
and  Holt's  model. 

A  slightly  modified  version  of  Kooker  and  Anderson's  model  is  also 
tested.  This  version  ignores  the  multiplier  Involving  internal  energy  on  the 
right  hand  side  of  Eq.  (3.14).  That  is,  for  this  case, 

P1  -1 

Teff{“*ec)  =  To  +  (TrT0)(1  *  5  *  t2^"c  i)i*  (3.16) 

This  run  is  made  to  see  the  effect  of  the  thermal  softening  term  on  the  com¬ 
puted  results.  It  is  also  made  since  we  are  comparing  these  results  to  those 
using  Carroll  and  Holt's  model,  which  does  not  include  such  a  term.  However, 
we  ran  the  code  up  to  t  =  35  ysec  with  no  detonation  occurring.  This  indeed 
demonstrates  that  with  the  thermal  softening  term,  the  resistance  to  pore  col¬ 
lapse  is  weaker. 


Table  6  Values  for  Parameters  Used  in  Eq.  (3.14)  and  (3.15) 
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Symbo 1 


Value  Assigned 
for  Hivelite  1086-A  (161 

3.44  x  IQ5  dyne-sec/cm 

2.54  x  IQ1-0  erg/g 

4.1023  x  IQ9  erg/g 

0.0 

2.758  x  108  dyne/cm2 
0.03  (non-dimensional) 
4.0  (non-dimensional) 


Table  7  Comparison  of  l$Qy  and  t$gy  using  Two  Different  Dynamic  Pore 
Collapse  Models.  The  inputs  are  p*  =  2  GPa,  t*  =  10  ysec,  a 
1.1176,  ax  =  0.5  mm.  For  Carroll  and  Holt's  [8]  model,  a  poPe 


radius  of  20  urn  was  used. 

Carroll  and  Holt 
Kooker  and  Anderson 
Difference 


L50T  (mm) 

tSDT  (,iSec) 

20.5032 

12.4 

15.3984 

30.0 

-  25% 

+  142% 
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